# -*- coding: utf-8 -*-
"""
Created on Wed Jun 11 11:37:53 2014

Field performance: takeoff, landing

@author: Maxim
"""
from numpy import zeros, ones, linspace, sin, cos, arcsin, cumsum, tan
import flight_conditions
from perf_tools import FlightMechanics, BasicInput


class Takeoff:
    """
    Takeoff analysis class. Takeoff consists of 3 segments: ground roll, 
    rotation and climbout. All these segments are analyzed separately. 
    Result of analysis is distance, time and fuel burned for each stage of 
    takeoff.
    
    Parameters
    ----------
    aircraft : aircraft object
    altitude : float, m
    """
    def __init__(self,aircraft,altitude=0.0,CLto=0.7,CLmaxTo=None,
                 hObstacle=15.0, breakFrictionCoef=0.04):
        if CLmaxTo==None:
            #CLmaxTo = aircraft.CLmax3D[0]
            CLmaxTo = aircraft.get_CLmax3D(25)
        self.CLto      = CLto
        self.CLmaxTo   = CLmaxTo
        self.hObstacle = hObstacle
        self.bi        = BasicInput(aircraft,CLmaxTo,altitude)
        self.fm        = FlightMechanics(aircraft)
        self.ac        = aircraft
        self.alt       = altitude
        self.atm       = flight_conditions.ISAtmosphere(altitude)
        self.rho       = self.atm.density
        self.mu        = breakFrictionCoef
        self._get_data()
        self.distance   = zeros(3)
        self.time       = zeros(3)
        self.fuelWeight = zeros(3)

    def _get_data(self):
        hw = self.ac.wing.leadingEdge[2] - self.ac.landingGear.groundContactZ[1]
        b  = self.ac.wing.span
        hw16 = (16.0*hw/b)**2.0
        self.Ge     = hw16 / (1.0+hw16) # ground effect
        self.Vstall = self.fm.get_Vstall()
        self.Vlo    = 1.1 * self.Vstall
        self.Vtr    = 1.2 * self.Vstall
        self.Vflare = 1.15 * self.Vstall
        self.thrustAv = self.ac.propulsion.get_thrust_available(self.alt)
        self.fuelFlow = self.ac.get_fuel_flow(self.Vstall,self.alt,self.thrustAv)

    def analyze(self,dtRot=1.5,powerSet=100.0):
        self.distance[0], self.time[0], self.fuelWeight[0] = self._ground_roll(powerSet)
        self.distance[1], self.time[1], self.fuelWeight[1] = self._rotation(powerSet,dtRot)
        self.distance[2], self.time[2], self.fuelWeight[2] = self._climbout(powerSet,dtRot)
        self.totalDistance = sum(self.distance)
        self.totalTime = sum(self.time)
        self.totalFuelWeight = sum(self.fuelWeight)
        return self.totalDistance, self.totalTime, self.totalFuelWeight

    def _ground_roll(self,powerSet,nSeg=25):
        tol  = 1e-5
        err  = tol+1.0
        g    = self.bi.g
        CDto = self.bi.Cd0 + self.Ge*self.bi.k*self.CLto**2.0
        V    = linspace(0.0,self.Vlo,nSeg)
        dV   = V[1]-V[0]
        Q    = 0.5 * self.bi.density * V**2.0
        L    = Q*self.bi.refArea*self.CLto
        D    = Q*self.bi.refArea*CDto
        Tav  = self.thrustAv
        ff   = self.fuelFlow
        T    = zeros(nSeg) + Tav
        ff   = zeros(nSeg) + ff
        dt   = ones(nSeg)
        while err>tol:
            tLast = sum(dt)
            Wf    = g*dV*cumsum(ff*dt)
            W     = self.bi.wt*ones(nSeg)-Wf
            F     = T-D-self.mu*(W-L)
            dt    = (W/g)/F
            err   = abs(tLast-sum(dt))
        ds = V*dt
        t = cumsum(dt*dV)
        s = cumsum(ds*dV)
        Sgr = s[-1]
        tgr = t[-1]
        Wfgr = W[0]-W[-1]
        return Sgr, tgr, Wfgr

    def _rotation(self,powerSet,dtRot=1.5):
        Srot = self.Vlo*dtRot
        wf = self.bi.g * self.fuelFlow*dtRot
        return Srot, dtRot, wf

    def _climbout(self,powerSet,dtRot=1.5):
        Wtr    = self.bi.wt - sum(self.fuelWeight)
        nFlare = 0.5*self.bi.density*self.Vflare**2.0 * self.bi.refArea*self.CLmaxTo / Wtr
        R      = self.Vflare / nFlare
        fuelFlow = self.ac.get_fuel_flow(self.Vflare,self.alt,self.thrustAv)
        Q  = self.bi.density*self.Vtr*self.Vtr / 2.0
        CL = Wtr / (Q*self.bi.refArea)
        CD = self.bi.get_drag_coefficient(CL)
        D  = Q*self.bi.refArea*CD
        climbAngle = arcsin((self.thrustAv-D)/Wtr)
        #-->
        Sflare = R*climbAngle
        hFlare = R*(1.0-cos(climbAngle))
        tFlare = Sflare / self.Vflare
        WfuelFlare = fuelFlow*tFlare
        #-->
        hTR  = self.hObstacle - hFlare
        RCtr = self.Vtr*sin(climbAngle)
        ttr  = hTR/RCtr
        Str  = self.Vtr*cos(climbAngle)*ttr
        WfTR = self.bi.g*fuelFlow*ttr

        Sclimb    = Sflare + Str
        timeClimb = tFlare + ttr
        WfClimb   = WfuelFlare + WfTR
        return Sclimb, timeClimb, WfClimb


# --- LANDING ---

class Landing():
    def __init__(self,aircraft,altitude,CLgr=0.7,CLmaxLDG=None,
                 hObstacle=15.0,breakFrictionCoef=0.40,reactionTime=1.0,loadFactor=1.2):
        if CLmaxLDG==None:
            #CLmaxLDG = aircraft.CLmax3D[1]
            CLmaxLDG = aircraft.get_CLmax3D(15)
        self.CLmaxLdg = CLmaxLDG
        self.hObst    = hObstacle
        self.bi       = BasicInput(aircraft,CLmaxLDG,altitude)
        self.fm       = FlightMechanics(aircraft)
        self.ac       = aircraft
        self.alt      = altitude
        self.atm      = flight_conditions.ISAtmosphere(altitude)
        self.rho      = self.atm.density
        self.mu       = breakFrictionCoef
        self.reactionTime = reactionTime

        self.idleThrust = 0.0
        self.loadFactor = loadFactor
        self.CLgr       = CLgr
        self.CLmaxLdg   = CLmaxLDG
        self._get_data()
        self.distance   = zeros(2)
        self.time       = zeros(2)
        self.fuelWeight = zeros(2)

    def _get_data(self):
        hw = self.ac.wing.leadingEdge[2] - self.ac.landingGear.groundContactZ[1]
        b  = self.ac.wing.span
        hw16 = (16.0*hw/b)**2.0
        self.Ge     = hw16 / (1.0+hw16) # ground effect
        self.Vstall = self.fm.get_Vstall()
        self.Vap = 1.3*self.Vstall
        self.Vtd = 1.1*self.Vstall
        self.Va  = 1.2*self.Vstall
        self.W   = self.bi.wt
        self.idleThrust = 0.0
        self.fuelFlow   = 0.0

    def analyze(self,nsegGR=25):
        self.distance[0],self.time[0],self.fuelWeight[0] = self._approach()
        self.distance[1],self.time[1],self.fuelWeight[1] = self._ground_roll(nsegGR)
        self.Wgr = self.bi.wt  -self.fuelWeight[0]*self.bi.g
        self.totalDistance   = sum(self.distance)
        self.totalTime       = sum(self.time)
        self.totalFuelWeight = sum(self.fuelWeight)
        return self.totalFuelWeight, self.totalTime, self.totalDistance

    def _get_cl(self,V):
        return 2.0*self.W / (self.atm.density*V*V*self.ac.wing.area)

    def _approach(self):
        T  = self.idleThrust
        ff = self.fuelFlow
        CL = self._get_cl(self.Va)
        CD = self.bi.get_drag_coefficient(CL)
        theta = arcsin(1.0/(CL/CD)-T/self.W)
        R     = self.Va**2.0 / (self.bi.g*(self.loadFactor - 1.0))
        hf = R*(1.0 - cos(theta))
        Sa = (self.hObst - hf) / tan(theta) + R*sin(theta)
        Vg = self.Va * cos(theta)
        ta = Sa / Vg
        Mfa = ff * ta
        return Sa, ta, Mfa*self.bi.g

    def _ground_roll(self,Nseg=25):
        tol = 1e-5
        err = tol+1.0
        g = self.bi.g
        CDldg = self.bi.Cd0 + self.Ge*self.bi.k * self.CLgr**2.0
        V     = linspace(self.Vtd,0.0,Nseg)
        dV    = V[1]-V[0]
        Q     = 0.5*self.bi.atm.density*V*V
        L     = Q*self.bi.refArea * self.CLgr
        D     = Q*self.bi.refArea * CDldg
        dt    = ones(Nseg)
        ff    = ones(Nseg)*self.fuelFlow
        while err>tol:
            tLast = sum(dt)
            Wf = self.bi.g * cumsum(ff*dt)
            W  = self.bi.wt*ones(Nseg) - Wf
            F  = self.idleThrust - D - self.mu*(W-L) #ma
            a  = F * g/ W
            dt = dV / a
            err = abs(tLast-sum(dt))
        ds = V*dt
        t  = cumsum(dt)
        s  = cumsum(ds)
        Sgr = s[-1] + self.reactionTime*self.Vtd
        tgr = t[-1] + self.reactionTime
        Wfgr = Wf[-1]
        Mfgr = Wfgr / self.bi.g
        return Sgr,tgr,Mfgr*g


def run_takeoff(ac,fieldAltitude=0.0,nsegGR=25):
    to = Takeoff(ac,fieldAltitude)
    data = to.analyze(nsegGR)
    return data[2],data[1],data[0] #fuel, time, distance

def run_landing(ac,fieldAltitude=0.0):
    ldg = Landing(ac,fieldAltitude)
    return ldg.analyze()


def field_performance_analysis(ac,detailedReport=False):
    fieldAlt = 0.0
    WfTO, tTO, RTO = run_takeoff(ac,fieldAlt)
    WfLDG, tLDG, RLDG = run_landing(ac,fieldAlt)
    if detailedReport:
        out = ''
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Takeoff distance',RTO,'m')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Takeoff fuel burned',WfTO,'kg')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Takeoff time',tTO,'sec')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Landing distance',RLDG,'m')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Landing fuel burned',WfLDG,'kg')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Landing time',tLDG,'sec')
    else:
        out = ''
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Takeoff distance',RTO,'m')
        out += '{:<23} = {:<12.1f}{:<6}\n'.format('Landing distance',RLDG,'m')
    return out
    

def run_test1():
    import aircraft_FW
    ac = aircraft_FW.load('Baseline1')
    print ac.get_mission_fuel()
    to = Takeoff(ac,0)
    print to.analyze()

    ldg = Landing(ac,0)
    print ldg.analyze()


if __name__=="__main__":
    run_test1()